function  output = processemg( EMGRAW, fs, LOWPASSRATE,NUMPASSES )
%%this function is created to process the emg data,the output should be the
%linear enveloped data 
%the input data should be raw emg data,samping rate,lowpass rate and the
%number of passes.

%% 读取数据
t = 0:1/fs:length(EMGRAW)*1/fs-1/fs; %这一步是根据采样频率计算要画出x轴的点的个数，也就是时间单位为秒
%画出原始EMG数据
figure();
subplot(5,1,1);
plot(t,EMGRAW);      %plot raw EMG data
xlabel('Time (s)');
%title('Raw EMG Data');
title('原始EMG数据');

%去除直流电 remove DC offset.
EMGDC = EMGRAW - mean(EMGRAW); %均值应该为0
subplot(5,1,2);
plot(t,EMGDC);
xlabel('Time (s)');
%title('EMG with DC Offset Removed');
title('去除直流电的EMG信号');

%滤波 filter the noise
Wn1 = 20/(fs/2); %the band of cutoff frequency
Wn2 = 450/(fs/2); %the band of cutoff frequency
[b,a]= butter(NUMPASSES,[Wn1 Wn2],'bandpass'); % 20-450 的带通滤波参数
EMGFILT = filtfilt(b,a,EMGDC); %带通滤波（去相位差后）
subplot(5,1,3);
plot(t,EMGFILT);
xlabel('Time (s)');
%title('EMG Signal with Noise Removed');
title('去除噪声的EMG信号');

%全波整流 full-wave rectification
EMGFWR = abs(EMGFILT); %全波整流，把正弦波的负值部分翻上去。
subplot(5,1,4);
plot(t,EMGFWR);
xlabel('Time (s)');
%title('EMG with Full-Wave Rectification');
title('全波整流EMG信号');

%线性包络 Linear Envelope
Wn = LOWPASSRATE/(fs/2);
[b,a] = butter(NUMPASSES,Wn,'low');
EMGLE = filtfilt(b,a,EMGFWR);% 低通滤波
subplot(5,1,5);
plot(t,EMGLE);
xlabel('Time (s)');
title('EMG信号的线性包络图');

output = EMGLE; % 输出结果
hold on
end